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Void coalescence and interaction in dynamic fracture of ductile metals have been investigated using 
three-dimensional strain-controlled multi-million atom molecular dynamics simulations of copper. 
The correlated growth of two voids during the coalescence process leading to fracture is investigated, 
both in terms of its onset and the ensuing dynamical interactions. Void interactions are quantified 
through the rate of reduction of the distance between the voids, through the correlated directional 
growth of the voids, and through correlated shape evolution of the voids. The critical inter-void lig- 
ament distance marking the onset of coalescence is shown to be approximately one void radius based 
on the quantification measurements used, independent of the initial separation distance between the 
voids and the strain-rate of the expansion of the system. The interaction of the voids is not reflected 
in the volumetric asymptotic growth rate of the voids, as demonstrated here. Finally, the practice 
of using a single void and periodic boundary conditions to study coalescence is examined critically 
and shown to produce results markedly different than the coalescence of a pair of isolated voids. 

PACS numbers: 61.72.Qq, 62.20.Mk, 62.20.Fe, 61.72.Lk 



I. INTRODUCTION 

The fracture of ductile metals at high strain rates has 
been understood at the microscopic level as a process of 
nucleation, growth and coalescence of voidsi^ Initially 
voids nucleate at the weak points in the material such 
as inclusions and/or grain boundary junctions. Once nu- 
cleated, the voids grow under the tensile stress, driven 
by the reduction in elastic energy. Eventually, the voids 
grow sufficiently large that they interact with each other, 
in some cases link through localized shear, coalesce into 
larger voids, and finally form the fracture surface>S»4 Con- 
siderable experimental and theoretical work has gone into 
the development of our understanding of fracture. The 
fracture process has been modeled at various levels, but 
most of the work involving the simulation of the activity 
of individual voids has concentrated on the void growth 
and its relationship to the plastic deformation of the sur- 
rounding material, including effects such as the localiza- 
tion of this shear deformation. Relatively little work has 
gone into the explicit modeling of void coalescence, and 
both the understanding of the physics of this process and 
the knowledge of how it should be implemented robustly 
in continuum fracture codes remain open issues. The 
point at which voids begin to coalesce during dynamic 
fracture is of considerable interest because complete frac- 
ture of the material typically ensues rapidly thereafter. 
As new experimental techniques have constrained the 
void growth models ever more stringently,^ a real need 
for a well developed theory of coalescence has arisen. 

Computationally void growth has been studied ex- 
tensively at the continuum level^SiLSiS also in dynamic, 
high strain-rate, conditions liSiii Recently we have stud- 
ied void growth at the atomistic level under high 
strain-rate expansion, motivated by spallation exper- 
iments The atomistic studies demonstrate 
that voids grow by emitting dislocations that carry away 



the material, platelets of atoms, from the void and are 
responsible for the plastic deformations needed to accom- 
modate significant void growth. There are also many 
recent studies of fracture in ductile metals with several 
holes or voids li^^iSiiSiSSiSi While these studies model the 
void growth explicitly with fairly sophisticated models of 
plasticity in many cases, they typically simplify the coa- 
lescence process to instantaneous unification of the voids 
when some threshold is reached, such as once the voids 
grow to within one diameter of each other. There are also 
studies that investigate the competition between void 
by void growth versus multiple void interaction in crack 
propagationiS^ Several earlier continuum studiesSiSiS^ 
and the one atomistic study known to us^fe of the coales- 
cence process have been typically conducted in effectively 
two dimensional and highly symmetric systems. 

This Article covers in detail a study of the onset of void 
coalescence. The first results of this study have been pre- 
sented in a LetteriSl Here we provide a more complete 
presentation of the results of our study of void coales- 
cence. In addition to a more detailed description of the 
results presented in the Letter, we describe different mea- 
sures of void interactions such as shape changes induced 
by a neighboring void, additional analysis of the behavior 
of the system including stress-strain curves and void vol- 
ume curves, and a new analysis of how the coalescence of 
isolated pairs of voids differs from that of voids in highly 
symmetric periodic arrays. 

In particular the goal of this Article is to quantify the 
point at which coalescence begins, as measured by a crit- 
ical inter-void ligament distance (ILD), and examine the 
mechanisms involved in the transition from independent 
void growth to coalescence. There are several ways in 
which two voids can interact. In the case of pure im- 
pingement, the voids only interact when they grow to 
the point that they intersect and join into a single void. 
In reality, the voids interact before they intersect. Their 



range of interaction is extended due to their elastic and 
plastic fields. Each void generates an elastic strain field 
of the form generally associated with centers of dilata- 
tion.^* The shear stress decreases with the distance from 
the void like r~^. For voids sufficiently close each void's 
growth rate is altered by the stress field of the proximal 
void. The modification of the elastic field can affect the 
initiation of plasticity, as well as the subsequent develop- 
ment of the plastic zone around the voids. The voids may 
interact through their plastic fields, too, in which case the 
fields may give rise to an increased hardening rate in a 
localized region or to thermal softening and shear local- 
ization. An argument due to Brown and Embury for a 
transition to shear deformation based on simple geomet- 
rical considerations suggests that the critical inter-void 
ligament distance, ILDc, should be equal to one diam- 
eter of a voidp^ that is, when the surfaces of a pair of 
voids are separated by one void diameter, they transi- 
tion from independent void growth to coalescence. It is 
at this point, they argue, that the dominant void pro- 
cess switches from the radial plastic flow around isolated 
growing voids to a large-scale shear deformation allowing 
the rapid coalescence of the pair of voids. However, more 
recent two-dimensional studies suggest that for distances 
between voids as large as six diameters the void growth 
rate is enhanced 

The use of atomistic techniques permits the analysis 
of the contributions of these competing mechanisms to 
the onset of void coalescence, as we describe in this Ar- 
ticle. We demonstrate the existence of, and compute, 
the critical inter-void ligament distance ILDc by starting 
with two voids well separated from each other and detect- 
ing the point at which correlated growth begins, marked 
both by the accelerated rate at which the two void sur- 
faces approach each other and by biased growth causing 
the voids to start to extend toward each other. These 
changes give an indication of the onset of the coalescence 
process, and it tests the argument by Brown and Em- 
burypS We also test the setup by Horstemeyer et alm& 
by varying the initial distances between the voids and 
measuring the asymptotic growth rate of the voids. The 
initial void-to- void distance below which the growth-rate 
is enhanced should give another candidate for the critical 
distance and measure it in a volumetric sense. It should 
be noted, however, that the three-dimensional void coa- 
lescence studied here with molecular dynamics, as indeed 
any 3D coalescence of roughly spherical voids, does not 
admit the two-dimensional shear mechanism proposed by 
Brown and Embury in its simplest form, and therefore 
these different analyses are not fully comparable. Also 
Horstemeyer et al. used a more symmetric setup than 
the simulations covered in this Article. 

The Article is organized as follows. The simulation 
method and the performed computations are introduced 
in Section ^ The basic mechanism, dislocation driven 
void growth, is demonstrated and the key reference pa- 
rameter, mean linear void size, is introduced in Sec- 
tion mil The interaction between voids are studied in 




FIG. 1: A sketch of the simulation configuration. The simula- 
tion box includes 120 x 120 x 120 FCC unit cells with periodic 
boundary conditions for a total of 6 912 000 atoms initially. 
From this two equal size voids are created by removing ap- 
proximately 3620 atoms for each void. The thin arrow shows 
the inter- void ligament distance (ILD), the shortest surface- 
to-surface distance between the voids, identified as A and B. 
The initial ILD ranges between one-half and five times the ini- 
tial void radius. The bold arrows denote triaxial, hydrostatic, 
expansion of the box. The strain-controlled expansion is ap- 
plied with constant strain-rates of e = 10** /sec and 10® /sec. 



TABLE I: The initial position for void B relative to void 
A and the center-to-center distances for various initial inter- 
void (surface-to-surface) ligament distances ILDo used in this 



study. 


ILDo 


position of void B 


separation between the 




relative to void A 


centers of voids A and B 


0.50 


[0.1334, 0.0622, 0.0290]L 


0.150L 


= 6.510 nm 


1.00 


[0.1778, 0.0829, 0.0387]L 


0.200L 


= 8.680 nm 


1.20 


[0.1956, 0.0912, 0.0426]L 


0.220L 


= 9.548 nm 


1.50 


[0.2223, 0.1037, 0.0483]L 


0.250L 


= 10.85 nm 


1.81 


[0.2500, 0.1166, 0.0544]L 


0.281L 


= 12.20 nm 


4.62 


[0.5000, 0.2332, 0.1087]L 


0.562L 


= 24.40 nm 



Section llVl Section llV Al uses two different distance mea- 
surements to study the interaction; Section IIV Bl intro- 
duces a shape parameter for the purpose. Finally the 
volume effects of the void growth are studied in Sec- 
tion El two separate voids are studied in Section IV Al 
and one void and the interaction with its periodic image 
are studied in Section [V Bl The paper ends with conclu- 
sions. Section IVll 
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II. METHOD AND SIMULATIONS 

We have performed a series of large-scale classical 
molecular dynamics (MD) simulations"^^ in single crys- 
tal face-centered cubic (FCC) systems using an empirical 
embedded-atom model (EAM) potential for copper i22iS 
The three dimensional (3D) simulation box consists of 
120 X 120 X 120 4-atom FCC cells with periodic bound- 
ary conditions for a total of 6 912 000 atoms, in most 
cases. In Section IV Bl where results from smaller system 
size simulations are presented, the details of the sizes are 
introduced. For the nearly 7 million atom simulations the 
MD code was parallelized using spatial domain decom- 
position and run in a massively parallel computer using 
from 64 to 256 processors. 

In the simulations the system is initially equilibrated 
using a thermostali^ at room temperature, T — 300 K, 
and a constant volume L'^ (with L — 43.3 nm) cho- 
sen to give ambient pressure, P ~ MPa. Once the 
system has reached equilibrium, atoms are removed to 
create two spherical voids in the system with radius 
rg = 0.05L = 2.17 nm: one centered in the box and 
the other located a distance ILDo-l-2ro away in the di- 
rection u = [0.8892054,0.41464327,0.19335135] from the 
first void. We refer to these as void A and void B, re- 
spectively, see Fig. ^ ILDq, the initial inter-void liga- 
ment distance, is the closest surface-to-surface distance 
between the voida^ and it is varied here, but the relative 
orientation of the voids is kept fixed. For ILDq we have 
used the values: 1.00, 1.20, 1.50, 1.81, and in some cases 
0.50 and 4.62. The unit for the ILD is the void diame- 
ter, d — 2r. The positions for the center of the void B 
and the distances between the voids are listed in Table 
Initially, the voids are equal in size, with approximately 
3620 atoms removed for each. This removal of atoms 
can be interpreted as an instantaneous debonding of two 
infinitely weakly bound inclusions. 

Once the voids are formed, the thermostat is turned 
ofi^, and dilatational strain is applied uniformly at a con- 
stant strain-rate e. The strain-controlled simulations^^ 
are carried out using the scaled coordinate formula- 
tion typically employed in the constant pressure method 
due to Parrinello and Rahman. Use of scaled coordi- 
nates prevents the spurious generation of elastic waves 
at the box boundaries. In this method the positions of 
the atoms are stored using rescaled coordinates between 
[0, 1). When calculating the forces and new positions of 
the atoms their coordinates are multiplied by a diagonal 
scaling matrix TL = {L^, Ly, L^}- This scaling matrix is 
updated each time step, when the expansion is applied, 
by multiplying the initial matrix TYq — {L, L, L} with 
the sum of the unit matrix and the strain matrix £ — t£, 
Hit) = 7-^0(1 + t^)- Applied strain-rates of e = 10^/sec 
and 10^/sec have been used with perfectly triaxial, or 
hydrostatic, expansion, £ — {Sx,Sy,ez} = {e,e,£}. Thus 
Lx{t) = Lyit) = Lz{t) = L{t) = ^1/3, where V is the 
volume of the box. A time step of 6.7 fs was used. More 
details of the simulation method can be found in Ref. 




FIG. 2: (Gray-scale) A three dimensional snapshot of the two 
voids with the prismatic dislocation loops forming from the 
voids. Only those atoms belonging to the void surfaces or to 
dislocation cores, stacking faults or other defects are shown 
(and a small number of extraneous atoms due to thermal fluc- 
tuations). The stacking fault ribbons are the broad plates be- 
tween leading and trailing partial dislocations. The snapshot 
is from the simulation with the initial inter-void ligament dis- 
tance ILDo=1.81 diameters, and the strain-rate e = 10^/sec. 
The strain at the snapshot is e = 2.93%. 

including analysis of growth of a single void of the same 
initial radius in non-triaxial expansion. 

III. DISLOCATIONS AND VOID GROWTH 

Let us start reviewing the simulation results by looking 
at some figures to visualize growth of the voids. While 
some void growth takes place through elastic stretching in 
the initial phases of the box expansion, significant void 
growth and void-void interaction take place only once 
plastic deformation has begun. The important role of 
plasticity leads us to consider in some detail how disloca- 
tions are generated and the effect of dislocation dynamics 
on void coalescence, which are the topics of this Section. 

In Fig. El a three-dimensional snapshot of the voids is 
shown from a simulation with the initial inter-void liga- 
ment distance ILDq =1.81 diameters and the strain-rate 
£ — 10^ /sec. In the plot only the atoms at crystallo- 
graphic defects such as void surfaces, dislocation cores, 
and stacking faults are shown. The decision of which 
atoms to plot is based on a geometrical criterion, a finite- 
temperature generalization of the centrosymmetry devi- 
ation.:i^ The power of this method of selecting atoms 
based on the centrosymmetry deviation is that the visu- 
alization can be done at finite temperature and on-the- 
fiy, so that the system need not to be cooled to zero- 
temperature and the atoms to be selected based on po- 
tential energy, where quenching can influence the sys- 
tem and prevent seeing the real configuration. The snap- 
shot in Fig. 12 is when the voids have started to grow by 
emitting dislocations, and thus the system has already 
evolved somewhat through plastic flow. The broad plates 
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FIG. 3: Dislocation activity at six instants of time shown on one particular slice through the system in order to expose the 
plasticity near the void surfaces. These snapshots are from the simulation with i = 10^ /sec and again only those atoms in 
dislocation cores, stacking faults, void surfaces or other defects are shown (see text for more details and Fig. |5|for a full three- 
dimensional figure). The dashed loop in panel (c) is drawn around a slice of a prismatic dislocation loop. The plane shown 
passes through the centers of both voids with normal [0.145, 0.145, -0.979]. The snapshots show the initial plasticity (a),(b), 
interacting plastic zones (c),(d) and the final coalescence (e),(f). The frames correspond to the following values of strain e: 
1.72%, 2.42%, 3.47%, 3.89%, 4.52%, and 5.21%, respectively. 



in the figure are stacking fault ribbons between leading 
and trailing partial dislocations. The first generation of 
dislocation loops have not yet totally formed and sepa- 
rated from the voids in this snapshot, but are about to 
do so, as can be seen as their first halves have already 
separated from the void surface. 

Figure 13 shows a series of visualizations of the crys- 
tal defects within a slice of width 4.5 A about a plane 
including centers of both voids at six different instants 
during coalescence. These snapshots are from the same 
simulation as in Fig. |21 and the same criterion for show- 
ing the atoms is used. Now only a slice is shown in order 
to reveal the dislocation activity near the void surface as 
the dislocation density increases. Figure |21 is a snapshot 
from a state of the system between panels Fig. Olb) and 
(c) . From the snapshots it is apparent that the deforma- 
tion mechanism involves the nucleation and propagation 
of dislocations, accommodating the void growth, and the 
interaction of the dislocationsi^ For example, the pris- 
matic dislocation loops punched out by the voids appear 
as roughly parallel line traces (due to the stacking fault 
ribbons) in the slice Fig. (SJc), as verified in the full 3D 
configuration (cf. Fig.[21l. Initially the dislocation activ- 



ity around each void is essentially symmetric [Figs. EJa) 
and (b)], as expected for independent void growth, but as 
the plastic fields evolve the void- void interaction is clearly 
evident both through interactions between the two plas- 
tic zones and bias due to the elastic fields [Fig. OIc)]. 
Once the dislocation density grows sufficiently high in 
the ligament region between the voids [Fig.|3{d)], void B 
begins to grow in the direction away from void A. Next 
the voids coalesce [Fig. Ole)], and continue to grow as 
one until ultimately the void coalesces with the its peri- 
odic images [subsequent to Fig. Elf)] so that the cavity 
percolates through the periodic system. The details of 
this final stage depend strongly on the periodic bound- 
ary conditions and will not be of interest here. 

The dislocation formation is closely related to the vol- 
ume evolution of the voids. The voids grow by emit- 
ting dislocation loops, driven by the reduction in the 
elastic energy as the increase in void volume allows the 
strained matrix material to relax. This relaxation can 
be detected from the saturation of the increasing stress 
in the system. All of these phenomena happen simul- 
taneously, see Ref. m In Fig. IJa) the mean (hy- 
drostatic) stress CTm is plotted with respect to strain e 
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FIG. 4: (a) Mean stress arn versus strain e from the simula- 
tions with the initial inter-void ligament distance ILDo=1.81 
and the strain-rates e = 10* /sec (dashed line) and 10® /sec 
(solid line), (b) Linear mean void size, f^^^ , (see the text for 
its definition) versus the strain from the same simulations as 
the data in (a). The void sizes are calculated until the co- 
alescence of the voids. The inset shows the stress from (a) 
versus linear mean void size from (b). Note the strain scales 
are different. 



(the control parameter in these simulations) for strain- 
rates e = lOVsec and lOVsec with ILDo==1.81. The 
stress is calculated with the virial expression Um — 

-W {^i + T.i,j\j>i • hj)-. where for atoms 

i and j the is the relative position, is the force, pi 
is the the momentum and is the mass (for a recent 
discussion of atomistic stress calculations see Ref. Isol) . 
Comparing Fig. E^a) with Fig. 01b) it is apparent that 
the start of the accelerated void growth due to plastic- 
ity is accompanied by, and indeed causes, the stress to 
plateau (at the same strain value) as the elastic dilata- 
tion is relieved. The void growth is shown in Fig.^b) by 
plotting linear void size f^^^ ^ where the single void frac- 
tion is / = Vvoid/V and V is the instantaneous volume of 



TABLE II: Critical strain Ec and critical linear mean void size 
fl^^ from Fig. |S] calculated as when the inter-void ligament 
distance crosses the line ILD=0.5d: fl'^ = /^/3(ILD=0.5d), 
= e(ILD=0.5d). 



ILDo 


e (sec ■"") 




ec (%) 


1.81 


10" 


0.153 


2.46 


1.81 


10® 


0.150 


3.48 


1.50 


10® 


0.132 


3.22 


1.20 


10® 


0.121 


3.06 


1.00 


10® 


0.106 


2.74 



the box at time t. The technique for calculating the void 
volume V^oid is described in Ref. The void growth, 
stress saturation and even the void coalescence take place 
at significantly smaller strains for slower strain-rates^ 
as can be seen from the figures, too. Therefore we con- 
clude that a natural way to plot quantities from different 
strain-rate simulations in the same figure is to plot them 
versus linear void size, /^/^. Plotting versus /^/'^ (a de- 
rived quantity) is preferred to plotting versus the strain 
because it reduces strain-rate effects. Thus by choosing 
ji/3 from Fig. Efb) as the reference quantity, the start 
of the deviation from the elastic behavior can be syn- 
chronized for the different strain-rates, see an example 
for the mean stress from the inset of Fig.^b), where the 
mean stress starts to saturate for both of the strain-rates 
at /^/^ ~ 0.09. In the two- void case we have chosen to 
use the mean void fraction /, which is calculated as the 
average void fraction of the two voids. 



Throughout the remainder of the Article we will typi- 
cally use the linear mean void size /^/^ as the reference 
quantity. The initial linear mean void fraction at ambient 

pressure in these simulations is J^^^ ~ (6.0 x lO"'')^/^ = 
0.084. The initial mean void size / = 6.0 x 10^'' is some- 
what larger than the value one gets from |7r(ro/L)^ = 
5.2 X lO"''. This potential source of confusion arises for 
two reasons: first, atoms that have their centers within 
the radius ro from the void center are removed for creat- 
ing the void. On the other hand, when the void volume 
is calculated the surface of the void is defined based on 
the centers of the remaining surface atoms. Second, the 
void surface relaxes somewhat after the void is formed. 
For more discussion of the void volume calculation, see 
Ref. For reference, the mean linear void size in Fig.|5] 
is /^/'^ = 0.111, and in the first four snapshots of Fig.|3 
/1/3 ^ 0.089^0.094, 0.149, and 0.195, respectively. After 
coalescence / is not measured. 
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linear mean void size (f ) 



FIG. 5: Evolution of the ILD and the critical ILD. (a) Dy- 
namical ILD, the distance between the surfaces of the voids 
along the line connecting the original center positions for vari- 
ous initial ILDq = 1.81, 1.50, 1.20, 1.00, plotted versus strain. 
For ILDo=1.81 the thick solid line and thick medium dashed 
line denote e — 10^/sec and 10*/sec, respectively. The thin 
solid lines show the hypothetical ILD for spherical voids with 
the same ILDo impinging freely on each other (see text). The 
short dashed line shows the hypothetical ILD computed by 
duplicating a single void (in the same box size as the two 
void simulations) at fixed centers (here the duplication is to 
the position with ILDo=1.81). (b) The same as in (a), but 
now plotted versus the linear average void size, and the 
distances are given in the units of the current average diam- 
eter d of the voids, calculated from their volumes assuming 
that they are spherical (see text). The horizontal line is at 
1LD=0.5 diameters, the value we identify as the ILDc. (Panel 
(b) is after Ref.0.) 



IV. THE INFLUENCE OF THE NEIGHBORING 
VOID 

A. Distance Measurements of the Voids 

Figure 13 offered several visual indications of the inter- 
action between voids. Clearly, the separation between 



the void surfaces (the ILD) serves as something akin to 
a reaction coordinate for the coalescence: the voids coa- 
lesce when it goes to zero. Other indications include the 
displacement of the center of a void as it grows prefer- 
entially toward the neighboring void and the change in 
the void growth rate as the voids interact. In this sec- 
tion we now quantify two of these elfects, the evolution 
of the ILD and the void center movement, in order to 
analyze the coalescence. In Section llV Bl we study the 
void shape evolution, and the void growth rate is studied 
m Section El 

In Figs.lSja) and (b) the dynamic evolution of the ILD 
has been plotted for strain-rates e = 10*/sec and 10^/sec 
and for various initial closest surface-to-surface distances 
between the voids ILDq. (In Figs. 12141 the data from the 
case ILDq ~ 1.81 were shown.) The dynamic ILD has 
been derived as the separation distance of the two sur- 
face atoms from the two voids that are closest to the 
line connecting the original centers of the voids. The 
raw data of the ILD have been plotted in Fig. Ela) and 
in Fig. Elb) scaled data are shown. For the horizontal 
axis of Fig.|S{b), the linear mean void size f-^^^ has been 
used in order to collapse data from different strain-rates 
£ — 10^ /sec and 10^ /sec in the same figure. For the 
vertical axis the inter-void ligament distance has been 
divided by the current diameter of a void, d (note, that 
it is not the initial value do = 2ro). The diameter d 
is calculated from the void volume assuming a spheri- 
cal shape (the formula is given below). This scaling is 
motivated by the ansatz that the void diameter sets the 
length scale for the system, and hence the relevant dis- 
tance between the voids is not the pure distance, but 
its ratio with the diameters of the voids. The void di- 
ameter affects both the elastic and plastic fields around 
the void. Initially the void separation distance decreases 
essentially smoothly until plasticity begins, eventually 
reaching zero. A transition occurs when the ILD starts 
to decrease noticeably faster than the free impingement 
line (the thin solid lines calculated from two indepen- 
dent spherical voids) indicating void interactions at the 
onset of coalescence. The thin solid lines are derived 
by calculating from the same simulations the average 
sizes of the two voids Vuoid and deriving the diameter 
as c? = 2 [^Vuoid) ■ Then the free impingement curve 
is given by ILD^^ee = (ILDq + do){l + te) — d. The ac- 
celerated ILD decrease associated with transition to coa- 
lescence takes place when the ILD reaches approximately 
one half: ILDc — 0.5 ±0.1 diameter or one radius, inde- 
pendently of ILDq or the strain-rate. Also note that a 
curve derived from a single void growth is provided to 
estimate the contribution of uncorrelated faceting effects 
(the "one void" curve at ILDo=1.81), and these effects 
are seen to be relatively small. They are most notice- 
able near the start of void growth where there is a rela- 
tively large upward fluctuation in the "one void" curve. 
The critical ILD of one radius is much lower than the 
Brown-Embury estimate, and it corresponds to a strain 
of 3.48% (/1/3 ~ 0.15) for ILDo=1.81 at e = lOVsec, 
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FIG. 6: Void center-of-mass displacement. The distance dcv 
from the original center of void to the instantaneous void cen- 
ter, projected onto the line connecting the original void cen- 
ters, is plotted versus the average void size to the point of 
coalescence (ILD ~ 0). The sign of the displacement dcv is 
positive for movement toward the other void. Solid and long 
dashed lines are for voids A and B, respectively. Strain-rates 
are i = 10^ /sec in (a) and e = 10* /sec in (b). The thin solid 
line is for a single void in the same size of the box and with 
the same radius and strain-rate projected to the same line. 
Here the distance dcv is given in the units of the original void 
diameter do- ILDo=1.81 in both (a) and (b). 



close to frame (c) of Fig. O The values for critical strain 
and linear mean void size when ILD=0.5d, derived from 
Fig. |3Jb) , are tabulated for the simulated systems in Ta- 
ble m In the very final stages the ligament is drawn 
under biaxial stress, and the flow switches from radial 
material transport to tangential transport as the mecha- 
nism switches from loop punching to drawing. This tran- 
sition is visible in Figs. [SJa) and (b) as slowing down of 
the reduction of ILD. At this point, the material is highly 
defective but it remains ductile. There is no abrupt frac- 
ture, as might be expected at larger length scales. Here 
the final coalescence involves an extended drawing and 
thinning of the ligament until rupture. 



Another measure of void interactions is whether the 
voids grow preferentially toward their neighbor. This ef- 
fect is quantified in Fig. |S1 which shows the movement 
of the center of mass of the void surface for the voids 
shown in Fig. |21(ILDo=1.81). We use here the center of 
mass of the void surface as the definition for calculations 
of the void center. The void surface has been derived by 
Voronoi triangulation based on the center points of the 
surface atoms (the same is done when the volume of the 
void is calculated) . See Ref . for the details of the void 
surface derivation. The distance 



dcv — 



,0 V 



1/2 



(1) 



is calculated between the original center of mass of the 
void surface [x^^, zj?^] and the projection [x' , y' , z'\ of 
the current center of void [xcd, y™, Zci^^ where the projec- 
tion is onto the line connecting the original void centers. 
The sign is positive if the void center has moved toward 
the other void and negative in the opposite case. 

Let us look first at the case with e = 10^/sec, Fig-EJa). 
After the void growth starts, the center of void A initially 
moves only slightly, but at about — 0.15 (ILD =0.5(i 
in Fig.[5Jb) and ]\^^ in Table ITU), it starts to move in the 
direction of the other void as the void growth becomes 
biased toward its neighbor. Just before coalescence the 
center of void A begins to move away from void B, as the 
growth is biased in the opposite direction. During this 
sequence, void B initially grows away from void A, then 
roughly in unison with void A (/^^^ = — 0.15) it be- 
gins to grow toward its neighbor, and before coalescence 
it too switches to growth away from the proximal void. 
This retrograde growth happens at the same point (after 

= 0.19) as the decrease of the ILD be gins to slow 
down in Fig.lSJb) [see also the snapshot in Fig.|3fd)]. The 
same phenomenon-first slow movement or repulsion from 
the void; then growth toward the nearby void at about 

= 0.15 (/y^)-holds in the e = lOVsec case, too. 
Fig. |SJb). However, the retrograde growth phenomenon 
is less pronounced in the e = 10^/sec case, as is the 
slowing down of the decrease of the ILD in Fig. [S{b). 
Indeed, only the growth toward the neighboring void is 
well above the noise in the e = 10®/sec case. As a refer- 
ence the movement of the center of a single void (in same 
box size) projected to the same line is plotted for both 
strain-rates, too. Comparing the single void case with 
the interacting voids with the same strain-rates, one sees 
that the maximum distance the centers of the interact- 
ing voids have moved is 2.5 to five times larger than the 
nanoscale random walk of the single void center, except 
for the void B in Fig. Efb). The noise is these curves 
appears to be dependent on the history, reflecting the 
nature of the plastic deformation processes involved. It 
is difficult to quantify the level of noise in such a non- 
Markov process, but it should be clear that the movement 
of the void center, especially the movement toward the 
neighboring void, is statistically significant and not just 
due to fluctuations at the void surface. 
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B. Shape Evolution of the Voids 

The presence of a nearby void not only affects the po- 
sition of the void but its shape, as well. The shape can be 
quantified by calculating multipole moments of spherical 
harmonics Q;™ ee ^ jYim{9, (/)) r"^ {9, (j)) dfl, see Refs.[ll 
and l4Ct Under fully triaxial expansion the void tends to 
an octahedral shape because of the index of the active 
glide planes in FCC crystals, the anisotropic elastic con- 
stants, and anisotropic surface energies. See also Ref. l4ll 
for a study of void shapes in FCC crystals. On the other 
hand, under uniaxial expansion the voids are of predomi- 
nantly ellipsoidal shapes aligned along the preferred axis, 
and they make a transition from a prolate to an oblate 
shape>i& References i42tt4Si44i have also considered oblate 
void shapes under uniaxial loading through continuum 
modeling. One may enquire whether the presence of a 
nearby void causes evolution to an ellipsoidal shape, see 
6.17. void A in Fig.l^tc). Ellipsoidal shapes can be quan- 
tified with quadrupole moments: 



^20 
Re Q21 
Im Q21 
Re Q22 
Im Q22 



(2) 



where — J r^{9,(p)dft is averaged over the surface 
of the void. In calculating the quadrupole moments the 
origin of the coordinates is taken to be the center of the 
void. In Eq. |2Jl the quadrupole moments are calculated 
with z as the preferred axis, whereas the physically pre- 
ferred axis is the line connecting the void centers. Thus, 
we transform the moments to the more natural coordi- 
nates using _D-matriceg^^ : 



Qi 
DL 



Z— — 



I rtvni 



{a,l3)Qin 



(3) 



where the Euler angle a defines the rotation between co- 
ordinates axes in (xy) plane (corresponding angle 4>) ; and 
the Euler angle /3 describes the rotation in z-axis (corre- 
sponding angle 9). dJ,j„/(/3) can be found from tables4& 
The quadrupole moments are plotted for the voids A and 
B from the simulation with ILDo=1.00 in Figs.[7|[a) and 
(b) , respectively. One sees from Fig. Cfa) that void A 
becomes markedly elliptical in the direction of the other 

void, as represented by Q20, when /^/^ ~ fc^^ ~ 0.106 
(from Table The quadrupole data from an iden- 

tical simulation but with only one void are plotted in 
Fig. [7fc) as a control, and a smaller variation in the 
quadrupole moment is observed simply due to fluctua- 
tions in the atomistic growth. Other cases with larger 
ILDq look about the same as Figs. [T^a) and (b) although 
the trend and especially the transition point may not 
be as clear. This variation may be due to the elastic 
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(b) void B 
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FIG. 7: (a) Quadrupole moments J^J of the surface of void 
A from the start of the simulation until coalescence. The 
coordinate axes have been rotated using Eq. ^ so that the 
z-axis is aligned with the center of void B. The initial inter- 
void ligament distance is ILDo=1.00 and the strain-rate is 
£ — 10^ /sec. (b) Quadrupole moments for the void B, lo- 
cated at [0.1778, 0.0829, 0.0387]L from void A. Now Eq. @ 
has been used so that the positive z-axs, is toward the void 
A. (c) Quadrupole moments (not rotated) for a single void 
in a simulation with same box size as in (a) and (b) and the 
strain-rate e = 10^/sec. 
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and plastic interactions causing random shape evolution 
longer before the critical ILD^, (ILD=0.5(i) is reached in 
cases when the voids are initially separated further from 
each other (larger ILDq). See Fig.[7Jc) for the single void 
case as an example of the random shape evolution. There 
— 0.09 is when the dislocation driven void growth 
starts, see inset of Fig.0Jb). Thus, in larger ILDq cases 
the random shape evolution suppresses the transition to 
the shape evolution due to the other void. 



V. VOID VOLUME EVOLUTION UNDER THE 
INFLUENCE OF THE SECOND VOID 

In Ref. 10, we briefly discussed how void growth after 
the onset of interaction but prior to coalescence is differ- 
ent from void growth for an isolated void. We address 
this further here. How does the correlated growth differ 
from the exponential growth of an isolated void?ii24i& 

A. A Pair of Voids 

We first examine the volume evolution for the same 
two- void simulations described above. Figurc|SJa) shows 
the void fraction, / = Vyoid/V , for the void A (also in the 
case of a single void in the same box size) with respect 
to the strain before coalescence. As can be seen from the 
figure the void grows as exp(200e), at least for the larger 
ILDo's and for the single void in the box. In Fig. IH^b) 
we have factored out the asymptotic growth rate from / 
in order to emphasize the differences between the curves 
and plotted versus linear void size f^^'^ . The void growth 
data for ILDo=4.62 and 1.81 coincide with the single void 
curve. The void growth rate with smaller ILDq's reach 
their asymptotic growth rate earlier. In the figure we 
have drawn as circles the void size values, where the dy- 
namic ILD's cross the line ILD=0.5(i in Fig. EJb). As 
can be seen from figure, there is no significant change in 
the void volume behavior when the voids start to inter- 
act. Therefore we conclude that the void growth rate is 
not affected by the interaction between the voids: thus 
the interaction cannot be detected through the growth 
rate. The key factor for the void growth rate is the rate 
at which the dislocations separate from the void, and it 
appears to remain unchanged in the vicinity of the sec- 
ond void. However, the location at the void surface from 
where the dislocation loops separate is affected by the in- 
teraction with another void, as seen in Section llVI in the 
accelerated reduction of ILD, in the movement of void 
center and in the shapes of the voids. 



B. A Single Void Interacting with its Periodic 
Images 

We have also performed a series of simulations of a sin- 
gle void with fixed initial radius size tq in various initial 




1.50 2.00 2.50 3.00 3.50 4.00 

Strain (e) [ % ] 



o 




FIG. 8: (a) Growth of the void A until coalescence presented 
by void fraction / = Vvoid/V versus strain for ILDo=0.50, 
1.00,1.20,1.50, 1.81, and 4.62 diameters as well as in a single 
void case (in the same box size) at e = 10^/sec. The asymp- 
totic behavior (before finite size effects) is exponential growth 
with exp(200e) as seen from the line drawn as a guide to the 
eye. (b) The void size / presented in (a) has been scaled with 
the exponential and plotted versus linear void size 

The 

circles point where the dynamical ILD's cross the horizontal 
line ILD=0.5d in Fig.l^b). 



box sizes Vq = I/^ in order to find the coalescence process 
of the void with its (six) periodic image(s), similar to the 
manner in which some continuum calculations of coales- 
cence have been done (cf. Refs.loland lisl) . The details of 
the different box sizes are reported in Table ITTll Figure |5| 
shows the data from this series of simulations. There we 
have scaled the void fraction / = V^oid/V with the expo- 
nential, exp(2G0e) (as in Fig.lSj) as well as {ro/L}^ in order 
to take into account different initial volumes Vq ~ L'^ oi 
the box. Also (rg / L) scaling is applied to the linear void 
f^^^ (horizontal axis). From the figure one concludes that 
the behavior is opposite to the results from the case of an 
isolated pair of voids, Section IVAI in an important way. 



10 



TABLE III: Sizes for simulated systems with various ro/L 
(ro = 2.17 nm is fixed) as the number of FCC cells, the equi- 
librium side length of the cube L at ambient pressure and 
room temperature, and the number of atoms in the box after 
the void is formed. The shortest initial inter-void ligament 
distance ILDo of the void with its periodic image in the units 
of the void diameter d is reported in the last column. 





FCC cells 


L 


atoms 


ILDo 


1/3 


18 X 18 X 18 


6.50 nm 


23328 


0.50 


2/9 


27 X 27 X 27 


9.75 nm 


78732 


1.25 


1/6 


36 X 36 X 36 


13.0 nm 


186624 


2.00 


1/8 


48 X 48 X 48 


17.3 nm 


442368 


3.00 


1/10 


60 X 60 X 60 


21.7 nm 


860396 


4.50 


1/20 


120 X 120 X 120 


43.3 nm 


6908379 


9.00 



a single void forces interstitial loops away from the void. 
This field decreases with the distance from the void as 
\/r^. At the near side of a neighboring void, this field 
would tend to drive interstitial loops into that void; in 
combination, it reduces that void's own stress field. At 
the far side, it adds to the other void's field, but the effect 
is smaller due to the greater distance. For a symmetric 
array, the effect is to reduce the maximum resolved shear 
stress across the surface and delay dislocation emission. 



VI. CONCLUSIONS 
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FIG. 9: Growth of a single void with varying initial box 
size. The initial radius of the void is kept constant ro = 
2.17nm, while the initial side length of the cube is varied as 
L =6.50nm, 9.75nm, 13.0nm, 17.3nm, 21.7nm, and 43.3nm. 
The void fraction / has been divided by (ro/L)^ in order to 
take into account different initial volumes Vo = oi the box. 
The vertical axis has been divided also with the exponential 
of the strain in order to show the asymptotic behavior as in 
Fig. |H] for the largest box size. The strain-rate is for all the 
simulations e — 10^ /sec. 



The smaller the box size, and hence the smaller the ILDo, 
the later the void starts to grow. This is easy to under- 
stand since measured from afar the two separate voids 
act as one big void, and thus grow faster than a smaller 
void. In a single void with its periodic image that picture 
does not hold. For example, there is no distance at which 
the stress field approaches a single void stress field. The 
periodic image only restricts the growth of the void. At a 
more mechanistic level, the net effect of the array of voids 
is to reduce the resolved shear stress driving dislocation 
emission from the void surface, whereas a single nearby 
void enhances this resolved shear stress on some regions 
of the void surface. An intuitive way to understand this 
phenomenon is to consider that the shear stress field of 



To summarize, interaction and coalescence of two voids 
in copper under tension have been simulated in multi- 
million-atom MD simulations. The effects of interac- 
tions between voids have been quantified by the increased 
reduction-rate of their separation, by the movement of 
their centers, and by their shape evolutions. The void 
interaction has also been visualized by detecting the dis- 
locations moving in the system using the generalized cen- 
trosymmetry parameter. The critical inter-void ligament 
distance has been found to be close to one void radius, in- 
dependent of the strain-rate or the initial separation dis- 
tance ILDq. The onset of coalescence occurs at the point 
that the plastic zones surrounding the voids first inter- 
act strongly. Signatures of coalescence have been found 
in the dynamic ILD curves and the void center move- 
ments, as explained in detail here, including reference to 
the stress-strain and void volume curves. A weaker sig- 
nature of the onset of coalescence has also been found 
in the void shape curves giving the quadrupole moment 
evolution. It has been demonstrated that the interaction 
of the voids is not reflected in the volumetric asymptotic 
growth rate of the voids. Finally, the coalescence process 
of an isolated pair of voids has been shown to be markedly 
different than the coalescence of a single void with its pe- 
riodic images, so the latter would not provide a reliable 
description of coalescence in typical low-symmetry con- 
figurations. 

In the future it would be interesting to study the cases 
with uniaxial expansion (and the various orientations of 
the voids with respect to the expansion direction) , differ- 
ent sizes of the voids relative to each other, other crystal 
structures as body-center cubic and hexagonal lattices, 
and systems including larger collections of voids. 
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